%Up biased transfn using E as reference fields
% Referece Sims et al 1971  Geophysics
% importance - the simple matrix notations l
% uses the matrix notations
%CHECKED WITH GEOTOOLS OK 30.4.3
% latest date 30.4.3

function[tf]=tf_up(SPMatrix),


if isstruct(SPMatrix) == 1, % for structur arryay
   nfreq = length(SPMatrix.spectra); % in momentory lapse of common sense i adopted structure
      for i = 1:nfreq,
   	SPM(i,:,:) = SPMatrix.spectra(i).data;
   end;
   SPM7=SPMatrix;
   clear SPMatrix;
   SPMatrix=SPM;
   
end;

A = size(SPMatrix);


if length(A) == 2,
   data(:,:)=SPMatrix;
   [EE,HE]=getmat(data);
   Z = inv(HE/EE); % inverse of admittance = impedance
   tf(1)=Z(1,2);% boring !! but i have to do this
   tf(2)=Z(2,1);
   tf(3)=Z(1,1);
   tf(4)=Z(2,2);
elseif length(A)==3,
   for i =1:A(1),
      data(:,:)=SPMatrix(i,:,:);
      [EE,HE]=getmat(data);
   Z = inv(HE/EE);
   tf(i,1)=Z(1,2);% boring !! but i have to do this
   tf(i,2)=Z(2,1);
   tf(i,3)=Z(1,1);
   tf(i,4)=Z(2,2);
end;
 
end;
tf=CompuErr(SPM7,tf); %variance need to be added 20.4.3
   
   
   
%helper function


function[EE,HE]=getmat(data),

j = sqrt(-1);

HxHx = data(1,1);
HyHy = data(2,2);
HzHz = data(3,3);
ExEx = data(4,4);
EyEy = data(5,5);

HxHy = data(2,1) - j*data(1,2);
HyHx = conj(HxHy);
HxEx = data(4,1) - j*data(1,4);
ExHx = conj(HxEx);
HxEy = data(5,1) - j*data(1,5);
EyHx = conj(HxEy);
ExHy = data(4,2) + j*data(2,4);
HyEx = conj(ExHy);
HyEy = data(5,2) - j*data(2,5);
EyHy = conj(HyEy);
ExEy = data(5,4)-j*data(4,5);
EyEx=conj(ExEy);

EE=[ExEx ExEy;EyEx EyEy];
HE=[HxEx HxEy;HyEx HyEy];

%---------------------------

function[tf]=CompuErr(SPM,tf),

nfreq=length(SPM.spectra);
for i = 1:nfreq,
   data(:,:) = SPM.spectra(i).data;
   c1 = cohe(tf(i,3),tf(i,1),data,'Ex');
   c2 = cohe(tf(i,2),tf(i,4),data,'Ey');
   dof = SPM.spectra(i).dof.avgt*SPM.spectra(i).dof.avgf;
   [tf(i,7),tf(i,5)] = BivError(c1,dof,data,'Ex');
   [tf(i,6),tf(i,8)] = BivError(c2,dof,data,'Ey');
end;
